Spectral Density Functionals for Electronic Structure Calculations 
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We introduce a functional of the local spectral electron density which can be used to to compute 
the total energy and the local spectral function of strongly-correlated materials. We illustrate the 
applicability of the method by using as an example the long-standing problem of the electronic 
structure of metallic plutonium. 
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Electronic structure approaches based on density func- 
tional theory (DFT) in its local density (LDA) and gen- 
eralized gradient (GGA) approximations |jj have suc- 
cessfully predicted the ground-state properties of weakly- 
correlated materials. Furthermore, the Kohn-Sham spec- 
trum of this approach has proved to be an excellent start- 
ing point for calculations of physical excitation energies 
using perturbative approaches such as the GW approxi- 
mation (2). 

On the other hand, this approach has failed qualita- 
tively in strongly-correlated systems such as Mott insula- 
tors, cuprates, manganites, and /-electron systems. Sev- 
eral directions have been pursued to treat such systems 
beyond LDA: The GW method ||, self- interaction cor- 
rected schemes ||, the LDA+U method Q and more 
recently LDA+DMFT The LDA+U approach 

combines the Dynamical Mean Field Theory (DMFT) 
which has been very successful in the context of 
model Hamiltonians, with state-of-the-art band struc- 
ture methods. This combination allows a coherent de- 
scription of the electron in real and momentum space 
and leads to the appearance of both Hubbard bands and 
quasiparticle bands. Hence, it is well suited to describe 
strongly-correlated materials. For recent applications of 
this method, see 

In this Letter we introduce a new functional formula- 
tion of the electronic structure problem which leads to the 
LDA+DMFT equations and we present an application of 
the full implementation of this method to plutonium. In 
our generalization of DMFT to the electronic structure 
problem, our basic idea is to introduce another relevant 
variable in addition to the density p, namely the local 
Green function. The latter is defined by projecting the 
full Green function onto a separate subset of correlated 
"heavy" orbitals distinguished by the index a from a com- 
plete set (indexed by a) of orbitals \a (r — R) = XaR of 
a tight-binding representation which we assume for sim- 
plicity to be orthogonal. The local Green function is 
therefore given by a matrix G with elements [jl~o|| 

G ab (iu>,R) = - (c aR (iuj)cf R (iuj)) = 
X * a (r - R) (</>(r, zc)V>+(r', iw)> Xb {r' - R)drdr'. (1) 



We then construct a functional T[p, G] which gives the 
exact free energy at a stationary point. This generalizes 
the effective action construction used in connection with 
either pure DMFT functionals jHJ or pure density func- 
tionals jnj . It amounts to first considering the partition 
function of the interacting electron gas in the presence 
of a static source coupled to the density and a dynamic 
source coupled to the local spectral function, and second, 
carrying out a Legendre transformation with respect to 
those sources. This functional can be constructed for- 
mally in perturbation theory, however its explicit form 
is not available just as in density functional theory. The 
success of DMFT in model hamiltonians suggests a use- 
ful approximation to the exact functional, namely the 
LDA+DMFT functional, which allows us to compute 
both total energies and spectra. 

To describe the new method, it is useful to introduce 
the notion of the Kohn-Sham potential Vks and its dy- 
namical analog E. They are defined as the functions that 
one needs to add to the kinetic energy matrix so as to 
obtain a given density and spectral function of the heavy 
orbitals namely: 



p(r)=T^2e l 



'(r|[i« + V 2 /2-Vi{-sr -E] _1 |r). (2) 



Vks is a function of r, the chemical potential p is set to 
zero throghout the paper, and E is given by 

E = E(r, r', iu) = £ x *(r - R)E a6 (^) Xb (r'-R), (3) 

abR 

Recall that we are truncating to the "heavy" subset {a} 
of orbitals. Then expression (0) becomes 



/>(*•) = T E E xkW h - ^ k - 



ib u.k 



X6k(r) 



(4) 



where H* b = (x a k| - V 2 + V K s\Xbk) is the one-electron 
hamiltonian in k-space. 
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The physical meaning of the dynamical potential £ is 
parallel to the meaning of the original Kohn-Sham po- 
tential Vks- h is the function that one needs to add to 
the correlated block of the one-electron hamiltonian in 
order to obtain the desired local Green function: 



E 

k 



[iw - H w - t{iu) 



ab 



(5) 



It is the frequency dependence of the dynamical potential 
which allows us to treat Hubbard bands and quasiparticle 
bands on the same footing. 

If the exact self energy of the problem is k- 
independent, then it coincides with X. The parallel ob- 
servation within DFT is: if the self energy of a model is 
k- and frequency- independent, then the self energy coin- 
cides with the Kohn-Sham potential. 

In terms of these quantities, the external potential and 
the matrix of local interactions U, we write down the 
DMFT+LDA functional: 

Flda+dmft(p> Vks,G, S) = 

—T e luja+ Tr log[iw + V 2 - V KS - S] 

- f %(r)p(r)(ir-^ e M0+ Tr[E( lW )GH] 

+ 1 V ext {v)p{r)dv+\ J ^^-dvdr' + E™ A [p] 



(0) 



<5[G] is the sum of the two-particle irreducible local di- 
agrams constructed with the local interaction matrix U, 
and the local heavy propagator G. $dc is the so-called 
double counting term which subtracts the average en- 
ergy of the heavy level already described by LDA. We 
write it here in its simple form, when only one Slater 
integral is included, i.e. <&dc = Un(n — l)/2 with 
n = T^.abGabiiuj)e luQ+ (see Ref. f§ for its general 
form). Expression [^ensures that the Greens function ob- 
tained from its extremization will satisfy the Luttinger 
theorem. 

An explicit form of 3>[G] is unavailable. We can calcu- 
late it by relating it to the free energy of an atom in a 
medium, W a u which we describe by an action S a t- To do 
this, we introduce the second central concept of DMFT: 
the Weiss field Gq 1 which represents the quadratic part 
of the action S a t- 



Sat[Go 1 ]= / ]T^(t)[G 1 (t,t')Wt') + 

ab 

J J2 U abcdci(T)c+(r)c c (T)c d (T). (7) 
T abed 



Here, G 1 = G at — A, where G a t describes the quadratic 
part of the action of the isolated atom, and A describes 



the bath surrounding the atom. The bath added to the 
atomic action produces the desired local Green function 
G a b- We can now relate the sum of local 



\CaC b )S a 

graphs $[(5] to the free energy of the atom in the medium 
W at = - log / exp[-S at ] via 



= W at [G- 



MGo 1 



G _1 )G — TrlogG. (8) 



The functional (g) can be viewed as a functional of four 
independent variables, since the stationary condition in 
the conjugate fields reproduces the definition of the dy- 
namical potential and the Weiss field. Extremizing it 
leads us to compute the Green function G a b(io;) Eq. (|5|) 
with the Kohn-Sham potential entering H^ b and with 



5<f> 



5$ 



DC 



SG a b(iuj) 8Gab{iw) 



(9) 



which identifies as the self energy of a generalized 

Anderson impurity model |?J in a bath. The "impurity" 
is characterized by a matrix of (local) levels 
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- H H ab 



(10) 



and is hybridized with the medium via the hybridization 
function A a &(z<x>). The matrix 5^Dc/SG a b(iuj) = 
U(n— 1/2) takes into account the double counting effects 
of the Coulomb interaction already contained in matrix 
H^ b . The hybridization function A a b(ioj) obeys the self- 
consistency condition 

iu - e ab - A ab (iw) - S afc (iw) = 



. k 



(ID 



The self-consistent DMFT loop with respect to A{iui) 
and Yi(iui) amounts to: (i) choosing some input hy- 
bridization A(iw), (ii) solving the Anderson impurity 
model with that hybridization and extracting the impu- 
rity Green function and self energy (iii) computing 
the on-site crystal Green function according to Eq. (j^), 
and (iv) recovering a new hybridization function A(iw) 
from Eq. ([□]). 

Since the charge density depends on E(zw), another 
self-consistent loop is required to obtain the extremum 
of Eq. (||) . This second loop updates the charge density 
once the DMFT loop delivers the self energy; this modi- 
fies the hopping integrals in the one-electron hamiltonian 
7? k . The new hoppings set up a new set of impurity lev- 
els e a b- We see that the solution of eqs.(p|)-(pTl) should be 
carried out in a double iterational loop: the DMFT loop 
finds T,(iui) for given p while the global one (DFT like) 
updates p. Within the DMFT technology, the implemen- 
tation of this global loop is an essential new development 
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of the present work which allows the determinatoion of a 
central quantity of the solid - its total energy. 

We now illustrate the method applied to the electronic 
structure of P u p| . Previous LDA HQ or LDA+U 
p5[ or model [fL6|] calculations have not been able to de- 
scribe the existence of both a compressed a phase and a 
low-density S phase, nor to provide a qualitative descrip- 
tion of its photoemission spectra (jl7|. We have recently 
reported the energy vs volume curve and the one-electron 
spectrum for a value of U which was estimated to cor- 
respond to that for S Pu Here, we study in detail 
the energy as a function of the parameter U and give 
a comparative analysis of the spectra in both a and 5 
phases We also report calculations on a bcc structure 
corresponding to the higher temperature e phase of Pu. 

We have implemented the spectral density functional 
approach described above using the linear muffin tin 
orbital (LMTO) method for electronic structure calcu- 
lations. We used a generalized non-orthogonal tight- 
binding LMTO representation jl8) which is most suited 
for describing the set of localized orbitals. The non- 
orthogonality corrections enter the definition for the 
Green function (||) in the form of an overlap matrix. 
To solve the impurity model for multi-orbital case of Pu 
we use a method which interpolates the self energy cal- 
culated from the Gutzwiller approximation Jl9| at low 
frequencies and a Hubbard- 1-like approximation |20[ at 
high frequencies |H| . Spin-orbit coupling effects are gen- 
erally important for actinide compounds and have been 
included in the calculation for Pu. The "full potential" 
terms have been neglected in the calculations through 
the utilization of the atomic sphere approximation with 
a one- kappa LMTO basis set p2fl . 

The necessary fc-space integrals for evaluating Green 
functions and charge densities have been carried out us- 
ing the tetrahedron method with on a reciprocal grid 
(8,8,8). The <5 and e phases of Pu have simple fee and bcc 
structures respectively, while the a phase has a compli- 
cated lattice with 16 atoms per unit cell. We determined 
the ground state energy as a function of the atomic vol- 
ume using fee and bcc geometries. Finally, the approach 
requires the effective Coulomb interaction between / elec- 
trons described by the Hubbard parameter U. Various 
estimates exist in the literature on the values of U in 
plutonium which indicate that the average interaction 
among / electrons is around 4 eV . 

We now report our calculated total energies for Pu in 
the fee structure. The temperature is fixed at 600-ftT, 
where the S phase is stable. To illustrate the importance 
of correlations we perform the calculations for several val- 
ues of U varying it from and 4 eV. The results are shown 
in Fig. 1: The C/=0 (GGA) curve shows a minimum at 
V/V$=0.72, close to the volume of the low-temperature 
a phase. Here Vg is the observed (5-phase volume. We 
expect that correlations should be less important for the 
compressed lattice in general, but there is no sign what- 



soever of the S phase in the U=0 calculation. Dramatic 
changes occur when U is about 4 eV. The curves show 
the possibility for a double minimum while the details de- 
pend sensitively on the actual value of U . We find that 
for U — 3.8 eV, the minimum occurs near V/V,5=0.80 
which corresponds to the volume of the a phase. When 
U increases by 0.2 eV the minimum occurs at V/V<5=1.05 
which corresponds to the volume of the 5 phase, in close 
agreement with experiment. Clearly, since the energies 
are so similar, we may expect that as temperature de- 
creases, the lattice undergoes a phase transition from the 
5 phase to the a phase with a remarkable decrease of the 
volume by 25%. 

We repeated our calculations for the bcc structure us- 
ing the temperature T = 900K where the e phase is sta- 
ble. Fig. 1 shows these results for U = 4 eV with a 
location of the minimum at around V/V,5=1.03. While 
the theory has a residual inaccuracy in determining the 
S and e phase volumes by a few percent, a hint of volume 
decrease with the S — > e transition is clearly reproduced. 
Thus, our first-principles calculations reproduce the main 
features of the experimental phase diagram of Pu. 

We now report our calculated spectral density of states 
for the fee structure using the volumes V/V,5=0.8 and 
V/V,5=1.05 corresponding to the a and 5 phases. Fig. 
2 shows the results of our dynamical mean-field calcula- 
tions as compared to the LDA-type calculations. We pre- 
dict the appearance of a strong quasiparticle peak near 
the Fermi level which exists both in the a and 5 phases. 
Recent advances have allowed the experimental determi- 
nation of these spectra, and our calculations are in accord 
with these measurements [0. The width of the quasi- 
particle peak in the a phase is found to be larger by 30 
per cent compared to the width in the S phase. This in- 
dicates that the low-temperature phase is more metallic, 
i.e. it has larger spectral weight in the quasiparticle peak 
and smaller weight in the Hubbard bands. 

In conclusion, we have developed a new method for 
electronic structure calculations which allows the simul- 
taneous evaluation of the total energy and the local elec- 
tronic spectral density. We applied the method to study 
the phase diagram of Pu and have reproduced several 
of its salient features. While we have not yet included 
the electronic and phonon entropy terms which are diffi- 
cult to evaluate numerically, many properties of realistic 
strongly correlated systems which are related to the total 
energy can now be explored using LDA+DMFT. 

The authors are indebted to E. Abrahams, A. Arko, 
J. Joyce and J. Thompson for helpful discussions. This 
research was supported by the US DOE, grant No. DE- 
FG02-99ER45761, by Los Alamos National Laboratory 
under subcontract No. 25124-001-012K, and by the NSF 
DMR-0096462. 
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FIGURE CAPTIONS 

Fig. 1. Total energy as a function of volume in Pu for 
different values of U calculated using the LDA+DMFT 
approach. Data for the fee lattice are computed at 
T=600K, while data for the bec lattice are given for 
T=900K. 

Fig. 2. Comparison between photoemission spec- 
tra of <5— Pu |l^] (circles) and calculated densities of 
states using the LDA+DMFT approach: the data for 
V/V s = 1.05, U = 4.0 eV (full line) and the data for 
V/Vs — 0.80, U — 3.8 eV (dashed line) correspond to the 
volumes of the S and a phases respectively. The result of 
the GGA calculation (dotted line) at V/Vs = 1 (U = 0) 
is also given. 
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Fig. 1 , Savrasov et. al 
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